Supplementary material
Area-level excess mortality in times of COVID-19 in Switzerland: geographical, socioeconomic and political determinants
Data
Age group | Sex | Observed | Expected (median) | Expected (lower bound) | Expected (upper bound) | Relative excess (median) | Relative excess (lower bound) | Relative excess (upper bound) |
|---|---|---|---|---|---|---|---|---|
40-59 | Female | 1,713 | 2,027 | 1,896 | 2,161 | 0.85 | 0.79 | 0.90 |
40-59 | Male | 2,966 | 2,542 | 2,390 | 2,699 | 1.17 | 1.10 | 1.24 |
60-69 | Female | 2,611 | 2,991 | 2,823 | 3,168 | 0.87 | 0.82 | 0.92 |
60-69 | Male | 4,478 | 3,675 | 3,478 | 3,911 | 1.22 | 1.14 | 1.29 |
70-79 | Female | 6,203 | 5,916 | 5,574 | 6,264 | 1.05 | 0.99 | 1.11 |
70-79 | Male | 8,972 | 6,901 | 6,534 | 7,276 | 1.30 | 1.23 | 1.37 |
80+ | Female | 27,541 | 20,790 | 19,817 | 21,719 | 1.32 | 1.27 | 1.39 |
80+ | Male | 20,292 | 20,274 | 19,263 | 21,436 | 1.00 | 0.95 | 1.05 |
Total | Total | 74,776 | 65,201 | 62,986 | 67,148 | 1.15 | 1.11 | 1.19 |
Models of observed and expected deaths by municipality
Step 1: iterative model development
To facilitate model development we start by only using the median excess mortality by municipality, age group and sex in 2020.
data1 = exp_deaths_2020_year %>%
group_by(canton, GMDNR, GMDNAME, age_group, id_space, sex, munici_observed, munici_pop,
density_high, density_low, across(starts_with("sep")), border, lang_fr, lang_it,
across(starts_with("vote")),type_urban,type_rural,type_periurban) %>%
summarise(munici_exp_deaths=median(munici_exp_deaths),
munici_excess=median(munici_excess)) %>%
mutate(E=log(ifelse(munici_exp_deaths==0,1e-4,munici_exp_deaths))) %>%
ungroup()
data2 = exp_deaths_2020_year %>%
group_by(canton, GMDNR, GMDNAME, age_group, id_space, sex, munici_observed, munici_pop,
density_high, density_low, across(starts_with("sep")), border, lang_fr, lang_it,
across(starts_with("vote")),type_urban,type_rural,type_periurban) %>%
summarise(munici_exp_deaths=mean(munici_exp_deaths),
munici_excess=mean(munici_excess)) %>%
mutate(E=log(ifelse(munici_exp_deaths==0,1e-4,munici_exp_deaths))) %>%
ungroup()
rm(exp_deaths_2020_year)
# gc()hyper.iid = list(theta = list(prior = "pc.prec", param = c(1, 0.01)))
hyper.bym2 = list(theta1 = list("PCprior", c(1, 0.01)),
theta2 = list("PCprior", c(0.5, 0.5)))
threads = parallel::detectCores()Model 1.0: no covariates
We use a model structure similar to Poisson regression, where \(O_{t,i,j,k}\), the number of observed deaths during week \(t\) in municipality \(i\), age group \(j\) and sex group \(k\), depends on the number of expected deaths \(E_{t,i,j,k}\) based on historical data and a log-linear predictor \(\log \lambda = \alpha + \beta X\).
\[ O_i \sim \text{Poisson}(\lambda E_i ) \\ \] At start, \(\lambda\) only includes one intercept parameter \(\alpha\), so that the estimate of \(\exp(\alpha)\) can be interpreted as an average relative excess mortality (that is, the ratio of observed on expected) for 2020. By adding covariates to \(\lambda\), we aim to disentangle the various factors that are associated with excess mortality at the local level.
We implement this model in R-INLA, a Bayesian inference package that
is especially adapted to spatial data. This is achieved in practice by
including \(\log (E_{i,j,k})\) as an
offset (although an alternative formulation based on the E
argument exists). During model development, we compare different model
versions based on the WAIC (lower values imply a better fit).
model1.0 = INLA::inla(munici_observed ~ 1 + offset(E),
data = data1,
family = "Poisson",
control.compute = list(config = TRUE, waic = TRUE),
quantiles = c(0.025, 0.5, 0.975),
num.threads = threads,
safe = TRUE)
summary(model1.0)Time used:
Pre = 2.2, Running = 0.645, Post = 0.121, Total = 2.96
Fixed effects:
mean sd 0.025quant 0.5quant 0.975quant mode kld
(Intercept) 0.161 0.004 0.154 0.161 0.168 0.161 0
Watanabe-Akaike information criterion (WAIC) ...: 38090.74
Effective number of parameters .................: 2.76
Marginal log-Likelihood: -19048.17
is computed
Posterior summaries for the linear predictor and the fitted values are computed
(Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)') mean 0.025quant 0.975quant
(Intercept) 1.175001 1.166609 1.183453[1] 1.175021As a sanity check, we find a relative excess mortality of 17.5% for 2020, that is coherent with a simple calculation (74,776 observed / 63,638 expected on average = 1.175). Remember that we excluded the age group 0-40, which explains why this is higher than numbers reported for Switzerland, generally around 10% for 2020. We can also look at the model fit and at the residuals. Obviously the model fit is not good here, as this basic model assumes a unique relative excess mortality for all areas, sexes and age groups.
Model 1.1: age and sex
We hypothesize that excess mortality affected different age and sex groups differently. We thus add the age group, the sex and the interaction of the two as covariates.
model1.1 = INLA::inla(munici_observed ~ - 1 + offset(E) +
sex:age_group,
data = data1,
family = "Poisson",
control.compute = list(config = TRUE, waic = TRUE),
quantiles = c(0.025, 0.5, 0.975),
num.threads = threads,
safe = TRUE)
summary(model1.1)Time used:
Pre = 1.37, Running = 0.762, Post = 0.146, Total = 2.28
Fixed effects:
mean sd 0.025quant 0.5quant 0.975quant mode kld
sexFemale:age_group40-59 -0.088 0.024 -0.136 -0.088 -0.041 -0.088 0
sexMale:age_group40-59 0.179 0.018 0.143 0.179 0.215 0.179 0
sexFemale:age_group60-69 -0.073 0.020 -0.111 -0.073 -0.034 -0.073 0
sexMale:age_group60-69 0.229 0.015 0.200 0.229 0.258 0.229 0
sexFemale:age_group70-79 0.063 0.013 0.038 0.063 0.087 0.063 0
sexMale:age_group70-79 0.300 0.011 0.280 0.300 0.321 0.300 0
sexFemale:age_group80+ 0.300 0.006 0.288 0.300 0.311 0.300 0
sexMale:age_group80+ 0.012 0.007 -0.001 0.012 0.026 0.012 0
Watanabe-Akaike information criterion (WAIC) ...: 36601.55
Effective number of parameters .................: 3.76
Marginal log-Likelihood: -18357.32
is computed
Posterior summaries for the linear predictor and the fitted values are computed
(Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)') mean 0.025quant 0.975quant
sexFemale:age_group40-59 0.9154638 0.8731221 0.9598588
sexMale:age_group40-59 1.1959522 1.1536771 1.2397765
sexFemale:age_group60-69 0.9299568 0.8949619 0.9663201
sexMale:age_group60-69 1.2575167 1.2212193 1.2948929
sexFemale:age_group70-79 1.0646097 1.0384432 1.0914355
sexMale:age_group70-79 1.3505051 1.3228475 1.3787409
sexFemale:age_group80+ 1.3494272 1.3335839 1.3654587
sexMale:age_group80+ 1.0124221 0.9985876 1.0264482[1] -1489.198As expected, the relative excess mortality varies a lot across age and sex groups. It’s very small in females aged 40-59 and 60-69 (in fact the data is compatible with no excess or even negative excess in both cases). It increases in females aged 70-79, and even more so aged 80+. It’s comparatively higher in males below 80, but somewhat surprisingly lower in males in age group 80+. This still corresponds to basic sanity checks with the data.
data1 %>%
group_by(sex,age_group) %>%
summarise(munici_observed=sum(munici_observed),
munici_exp_deaths=sum(munici_exp_deaths)) %>%
mutate(ratio=munici_observed/munici_exp_deaths)# A tibble: 8 × 5
# Groups: sex [2]
sex age_group munici_observed munici_exp_deaths ratio
<chr> <chr> <int> <dbl> <dbl>
1 Female 40-59 1713 1870. 0.916
2 Female 60-69 2611 2807 0.930
3 Female 70-79 6203 5826 1.06
4 Female 80+ 27541 20409 1.35
5 Male 40-59 2966 2480. 1.20
6 Male 60-69 4478 3560. 1.26
7 Male 70-79 8972 6643 1.35
8 Male 80+ 20292 20042. 1.01 We observe an improvement of the model fit, not easy to spot on the plot because of the large number of points, but made clear by the large decrease in WAIC.
Model 1.2: spatial variability
We now account for spatial variability, first in a simple way using an i.i.d. random effect, so that all municipalities can vary independently from each other around a global average. Note that this “municipality effect” applies the same to all age and sex groups.
model1.2 = INLA::inla(munici_observed ~ - 1 + offset(E) +
sex:age_group +
f(id_space, model = "iid"),
data = data1,
family = "Poisson",
control.compute = list(config = TRUE, waic = TRUE),
quantiles = c(0.025, 0.5, 0.975),
num.threads = threads,
safe = TRUE)
summary(model1.2)Time used:
Pre = 2.32, Running = 2.1, Post = 0.287, Total = 4.71
Fixed effects:
mean sd 0.025quant 0.5quant 0.975quant mode kld
sexFemale:age_group40-59 -0.086 0.024 -0.133 -0.086 -0.039 -0.086 0
sexMale:age_group40-59 0.182 0.018 0.145 0.182 0.218 0.182 0
sexFemale:age_group60-69 -0.070 0.020 -0.108 -0.070 -0.031 -0.070 0
sexMale:age_group60-69 0.231 0.015 0.202 0.231 0.261 0.231 0
sexFemale:age_group70-79 0.065 0.013 0.040 0.065 0.090 0.065 0
sexMale:age_group70-79 0.303 0.011 0.282 0.303 0.324 0.303 0
sexFemale:age_group80+ 0.303 0.006 0.291 0.303 0.315 0.303 0
sexMale:age_group80+ 0.015 0.007 0.001 0.015 0.029 0.015 0
Random effects:
Name Model
id_space IID model
Model hyperparameters:
mean sd 0.025quant 0.5quant 0.975quant mode
Precision for id_space 4628.86 3952.99 1559.50 3411.26 14357.23 2488.85
Watanabe-Akaike information criterion (WAIC) ...: 36576.63
Effective number of parameters .................: 12.66
Marginal log-Likelihood: -18353.31
is computed
Posterior summaries for the linear predictor and the fitted values are computed
(Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)') mean 0.025quant 0.975quant
sexFemale:age_group40-59 0.9175715 0.8750433 0.9621683
sexMale:age_group40-59 1.1991021 1.1565772 1.2431938
sexFemale:age_group60-69 0.9325095 0.8973085 0.9690939
sexMale:age_group60-69 1.2601251 1.2236144 1.2977290
sexFemale:age_group70-79 1.0672440 1.0408572 1.0943054
sexMale:age_group70-79 1.3535962 1.3256746 1.3821142
sexFemale:age_group80+ 1.3535281 1.3371410 1.3701555
sexMale:age_group80+ 1.0149133 1.0007999 1.0292418[1] -24.91934The age and sex effect remains similar, but the model fit as measured by the WAIC is improved now that we account for local differences.
We find noisy estimates in some places, suggesting issues related to small area estimation. One solution is to partially pool information between municipalities that are geographically linked based on a spatial structure.
Model 1.3: structured spatial variability
We still focus on spatial variability, but now the municipalities are no longer independent: we account for the correlation between neighboring municipalities with a BYM model. Neighboring municipalities are defined as municipalities sharing a border. This will allow us to differentiate between what can be attributed to a municipality in particular, and what can be attributed to regional effects (like a COVID wave).
model1.3 = INLA::inla(munici_observed ~ - 1 + offset(E) +
sex:age_group +
f(id_space, model = "bym2", graph = "data/nb/gg_wm_q.adj", scale.model = TRUE,
hyper = hyper.bym2, constr=TRUE),
data = data1,
family = "Poisson",
control.compute = list(config = TRUE, waic = TRUE),
quantiles = c(0.025, 0.5, 0.975),
num.threads = threads,
safe = TRUE)
summary(model1.3)Time used:
Pre = 20.1, Running = 5.28, Post = 0.503, Total = 25.9
Fixed effects:
mean sd 0.025quant 0.5quant 0.975quant mode kld
sexFemale:age_group40-59 -0.084 0.024 -0.131 -0.084 -0.036 -0.084 0
sexMale:age_group40-59 0.183 0.018 0.147 0.183 0.220 0.183 0
sexFemale:age_group60-69 -0.068 0.020 -0.106 -0.068 -0.029 -0.068 0
sexMale:age_group60-69 0.235 0.015 0.205 0.235 0.264 0.235 0
sexFemale:age_group70-79 0.068 0.013 0.042 0.068 0.093 0.068 0
sexMale:age_group70-79 0.306 0.011 0.285 0.306 0.327 0.306 0
sexFemale:age_group80+ 0.305 0.006 0.293 0.305 0.317 0.305 0
sexMale:age_group80+ 0.017 0.007 0.003 0.017 0.032 0.017 0
Random effects:
Name Model
id_space BYM2 model
Model hyperparameters:
mean sd 0.025quant 0.5quant 0.975quant mode
Precision for id_space 1010.491 243.148 623.029 980.115 1574.598 919.535
Phi for id_space 0.953 0.045 0.829 0.967 0.996 0.989
Watanabe-Akaike information criterion (WAIC) ...: 36435.70
Effective number of parameters .................: 10.33
Marginal log-Likelihood: -17455.00
is computed
Posterior summaries for the linear predictor and the fitted values are computed
(Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)') mean 0.025quant 0.975quant
sexFemale:age_group40-59 0.9195800 0.8768831 0.9643569
sexMale:age_group40-59 1.2012016 1.1585014 1.2454773
sexFemale:age_group60-69 0.9346884 0.8993344 0.9714335
sexMale:age_group60-69 1.2645985 1.2278356 1.3024645
sexFemale:age_group70-79 1.0698472 1.0432740 1.0971007
sexMale:age_group70-79 1.3579398 1.3297173 1.3867668
sexFemale:age_group80+ 1.3565543 1.3398154 1.3735241
sexMale:age_group80+ 1.0174590 1.0030824 1.0320515[1] -140.9299We see that the structure accounts for a large part of the spatial variability (Phi estimated to 95%). This addition also improves the model fit as measured by the WAIC. The following map allows to look at specific municipality effects.
We observe that many of the municipalities with higher relative excess mortality are in the western and southern parts, the ones that were hit first by COVID-19 in 2020. We also observe areas with higher excess in the North and Northeastern parts. These largely correspond to areas that were hit the most during the first and the second COVID-19 waves of spring and fall 2020 (Konstantinoudis et al. 2022).
Model 1.4: local characteristics
Having accounted for regional variability (arguably caused by COVID-19 waves of different timings and scales), we move on to explore the effect of local characteristics at the municipality level.
Rural/urban
The Federal Statistical Office classifies Swiss municipalities in 3 classes: urban, peri-urban or rural (https://www.bfs.admin.ch/bfs/en/home/statistics/territory-environment/nomenclatures/gemtyp.html). We add this covariate to the model taking the “rural” category as the reference.
model1.4a = INLA::inla(munici_observed ~ - 1 + offset(E) +
sex:age_group +
f(id_space, model = "bym2", graph = "data/nb/gg_wm_q.adj", scale.model = TRUE,
hyper = hyper.bym2, constr=TRUE) +
type_periurban + type_urban,
data = data1,
family = "Poisson",
control.compute = list(config = TRUE, waic = TRUE),
quantiles = c(0.025, 0.5, 0.975),
num.threads = threads,
safe = TRUE)
summary(model1.4a)Time used:
Pre = 21.9, Running = 5.18, Post = 0.533, Total = 27.6
Fixed effects:
mean sd 0.025quant 0.5quant 0.975quant mode kld
type_periurban -0.016 0.013 -0.041 -0.016 0.009 -0.016 0
type_urban -0.033 0.011 -0.054 -0.033 -0.011 -0.033 0
sexFemale:age_group40-59 -0.061 0.026 -0.112 -0.061 -0.011 -0.061 0
sexMale:age_group40-59 0.205 0.020 0.166 0.205 0.245 0.205 0
sexFemale:age_group60-69 -0.045 0.021 -0.087 -0.045 -0.002 -0.045 0
sexMale:age_group60-69 0.257 0.017 0.223 0.257 0.290 0.257 0
sexFemale:age_group70-79 0.090 0.015 0.060 0.090 0.120 0.090 0
sexMale:age_group70-79 0.328 0.014 0.302 0.328 0.355 0.328 0
sexFemale:age_group80+ 0.328 0.011 0.307 0.328 0.349 0.328 0
sexMale:age_group80+ 0.040 0.011 0.018 0.040 0.062 0.040 0
Random effects:
Name Model
id_space BYM2 model
Model hyperparameters:
mean sd 0.025quant 0.5quant 0.975quant mode
Precision for id_space 1066.028 259.774 653.141 1033.240 1669.611 968.02
Phi for id_space 0.954 0.044 0.834 0.968 0.996 0.99
Watanabe-Akaike information criterion (WAIC) ...: 36429.58
Effective number of parameters .................: 9.95
Marginal log-Likelihood: -17466.15
is computed
Posterior summaries for the linear predictor and the fitted values are computed
(Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)') mean 0.025quant 0.975quant
type_periurban 0.9838908 0.9597851 1.0086042
type_urban 0.9677412 0.9472574 0.9886839
sexFemale:age_group40-59 0.9404828 0.8942800 0.9890729
sexMale:age_group40-59 1.2280313 1.1801812 1.2778216
sexFemale:age_group60-69 0.9563903 0.9169689 0.9975065
sexMale:age_group60-69 1.2924242 1.2496664 1.3366456
sexFemale:age_group70-79 1.0944904 1.0619653 1.1280120
sexMale:age_group70-79 1.3886530 1.3519699 1.4263322
sexFemale:age_group80+ 1.3883682 1.3596361 1.4177091
sexMale:age_group80+ 1.0407590 1.0182590 1.0637570[1] -6.11736
On average, urban and to a lesser extent peri-urban municipalities appear to have a lower excess mortality than municipalities classified as rural.
Socio-economic position
The Swiss neighbourhood index of socio-economic position provides an estimate of socio-economic position (SEP) based on census data for 1.5 million buildings (Panczak et al. 2023). We consider the median index of each municipality, then group municipalities in quintiles before adding to the model (reference is 5th quintile with highest SEP).
model1.4b = INLA::inla(munici_observed ~ - 1 + offset(E) +
sex:age_group +
f(id_space, model = "bym2", graph = "data/nb/gg_wm_q.adj", scale.model = TRUE,
hyper = hyper.bym2, constr=TRUE) +
sep1 + sep2 + sep3 + sep4,
data = data1,
family = "Poisson",
control.compute = list(config = TRUE, waic = TRUE),
quantiles = c(0.025, 0.5, 0.975),
num.threads = threads,
safe = TRUE)
summary(model1.4b)Time used:
Pre = 20.7, Running = 5, Post = 0.506, Total = 26.3
Fixed effects:
mean sd 0.025quant 0.5quant 0.975quant mode kld
sep1 0.065 0.016 0.034 0.065 0.096 0.065 0
sep2 0.037 0.013 0.012 0.037 0.062 0.037 0
sep3 0.031 0.012 0.007 0.031 0.055 0.031 0
sep4 0.004 0.011 -0.018 0.004 0.025 0.004 0
sexFemale:age_group40-59 -0.106 0.025 -0.156 -0.106 -0.057 -0.106 0
sexMale:age_group40-59 0.160 0.020 0.121 0.160 0.199 0.160 0
sexFemale:age_group60-69 -0.090 0.021 -0.131 -0.090 -0.048 -0.090 0
sexMale:age_group60-69 0.212 0.017 0.179 0.212 0.245 0.212 0
sexFemale:age_group70-79 0.045 0.015 0.016 0.045 0.074 0.045 0
sexMale:age_group70-79 0.284 0.013 0.258 0.284 0.309 0.284 0
sexFemale:age_group80+ 0.283 0.010 0.264 0.283 0.302 0.283 0
sexMale:age_group80+ -0.005 0.010 -0.025 -0.005 0.015 -0.005 0
Random effects:
Name Model
id_space BYM2 model
Model hyperparameters:
mean sd 0.025quant 0.5quant 0.975quant mode
Precision for id_space 1391.728 391.616 791.805 1335.342 2320.703 1225.723
Phi for id_space 0.941 0.064 0.762 0.963 0.997 0.993
Watanabe-Akaike information criterion (WAIC) ...: 36427.72
Effective number of parameters .................: 9.39
Marginal log-Likelihood: -17474.90
is computed
Posterior summaries for the linear predictor and the fitted values are computed
(Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)') mean 0.025quant 0.975quant
sep1 1.0672165 1.0347536 1.1005102
sep2 1.0380694 1.0122114 1.0644357
sep3 1.0318276 1.0075048 1.0566084
sep4 1.0038915 0.9825334 1.0256664
sexFemale:age_group40-59 0.8990845 0.8554782 0.9449308
sexMale:age_group40-59 1.1736932 1.1288780 1.2203223
sexFemale:age_group60-69 0.9143268 0.8775121 0.9527100
sexMale:age_group60-69 1.2356631 1.1957205 1.2769916
sexFemale:age_group70-79 1.0459372 1.0162552 1.0765405
sexMale:age_group70-79 1.3279471 1.2948206 1.3620099
sexFemale:age_group80+ 1.3268964 1.3024788 1.3519439
sexMale:age_group80+ 0.9950145 0.9754383 1.0150904[1] -7.977769
A gradient appears clearly, so that we can conclude that municipalities of lowest median SEP had a higher relative excess mortality in 2020 compared to municipalities of highest median SEP (5th quintile).
International borders
We now consider whether the municipality belongs to a cross-border labor region (https://www.bfs.admin.ch/bfs/de/home/grundlagen/raumgliederungen.assetdetail.8706500.html). This identifies municipalities with a high level of connectivity with neighboring countries France, Germany and Italy.
model1.4c = INLA::inla(munici_observed ~ - 1 + offset(E) +
sex:age_group +
f(id_space, model = "bym2", graph = "data/nb/gg_wm_q.adj", scale.model = TRUE,
hyper = hyper.bym2, constr=TRUE) +
border,
data = data1,
family = "Poisson",
control.compute = list(config = TRUE, waic = TRUE),
quantiles = c(0.025, 0.5, 0.975),
num.threads = threads,
safe = TRUE)
summary(model1.4c)Time used:
Pre = 19.9, Running = 4.46, Post = 0.514, Total = 24.9
Fixed effects:
mean sd 0.025quant 0.5quant 0.975quant mode kld
border 0.025 0.014 -0.003 0.025 0.053 0.025 0
sexFemale:age_group40-59 -0.088 0.024 -0.136 -0.088 -0.041 -0.088 0
sexMale:age_group40-59 0.179 0.019 0.143 0.179 0.216 0.179 0
sexFemale:age_group60-69 -0.072 0.020 -0.111 -0.072 -0.033 -0.072 0
sexMale:age_group60-69 0.230 0.015 0.200 0.230 0.260 0.230 0
sexFemale:age_group70-79 0.063 0.013 0.037 0.063 0.089 0.063 0
sexMale:age_group70-79 0.301 0.011 0.280 0.301 0.323 0.301 0
sexFemale:age_group80+ 0.300 0.007 0.287 0.300 0.314 0.300 0
sexMale:age_group80+ 0.013 0.008 -0.002 0.013 0.028 0.013 0
Random effects:
Name Model
id_space BYM2 model
Model hyperparameters:
mean sd 0.025quant 0.5quant 0.975quant mode
Precision for id_space 1020.153 246.177 628.234 989.282 1591.621 927.702
Phi for id_space 0.951 0.047 0.823 0.966 0.996 0.989
Watanabe-Akaike information criterion (WAIC) ...: 36434.49
Effective number of parameters .................: 10.75
Marginal log-Likelihood: -17461.15
is computed
Posterior summaries for the linear predictor and the fitted values are computed
(Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)') mean 0.025quant 0.975quant
border 1.0252579 0.9973396 1.0540144
sexFemale:age_group40-59 0.9153587 0.8726247 0.9601875
sexMale:age_group40-59 1.1961389 1.1532780 1.2405960
sexFemale:age_group60-69 0.9305329 0.8950596 0.9674148
sexMale:age_group60-69 1.2590251 1.2219416 1.2972396
sexFemale:age_group70-79 1.0649608 1.0379928 1.0926363
sexMale:age_group70-79 1.3517208 1.3228403 1.3812440
sexFemale:age_group80+ 1.3503193 1.3323419 1.3685766
sexMale:age_group80+ 1.0127815 0.9976015 1.0282112[1] -1.206236
We observe higher relative excess mortality in municipalities in the vicinity of international borders, although the data remains compatible with no difference.
Language
Most Swiss municipalities have one official language: German, French or Italian. A few municipalities have several official languages, but given the relatively low numbers, we consider only the majority language. The difficulty is the colinearity between language regions and the first COVID-19 wave of 2020, that primarily affected Ticino (Italian) and Southwestern Switzerland (French), mostly because of how the initial global spread of COVID-19 occurred (with large early epidemics in Italy and then France). These effects are much larger than any effect that could be attributed to cultural differences between language regions, so it is quite difficult to estimate the latter. We still attempt to do so by adding the language of each municipality (reference is German) to our model.
model1.4d = INLA::inla(munici_observed ~ - 1 + offset(E) +
sex:age_group +
f(id_space, model = "bym2", graph = "data/nb/gg_wm_q.adj", scale.model = TRUE,
hyper = hyper.bym2, constr=TRUE) +
lang_fr + lang_it,
data = data1,
family = "Poisson",
control.compute = list(config = TRUE, waic = TRUE),
quantiles = c(0.025, 0.5, 0.975),
num.threads = threads,
safe = TRUE)
summary(model1.4d)Time used:
Pre = 20.1, Running = 5.66, Post = 0.511, Total = 26.3
Fixed effects:
mean sd 0.025quant 0.5quant 0.975quant mode kld
lang_fr 0.082 0.012 0.057 0.082 0.106 0.082 0
lang_it 0.135 0.021 0.094 0.136 0.175 0.136 0
sexFemale:age_group40-59 -0.115 0.024 -0.163 -0.115 -0.067 -0.115 0
sexMale:age_group40-59 0.152 0.019 0.115 0.152 0.189 0.152 0
sexFemale:age_group60-69 -0.100 0.020 -0.139 -0.100 -0.061 -0.100 0
sexMale:age_group60-69 0.204 0.015 0.174 0.204 0.234 0.204 0
sexFemale:age_group70-79 0.035 0.013 0.009 0.035 0.061 0.035 0
sexMale:age_group70-79 0.274 0.011 0.252 0.274 0.296 0.274 0
sexFemale:age_group80+ 0.272 0.007 0.258 0.272 0.287 0.272 0
sexMale:age_group80+ -0.015 0.008 -0.031 -0.015 0.001 -0.015 0
Random effects:
Name Model
id_space BYM2 model
Model hyperparameters:
mean sd 0.025quant 0.5quant 0.975quant mode
Precision for id_space 4667.906 2919.752 1458.418 3923.820 12385.132 2854.588
Phi for id_space 0.819 0.158 0.396 0.868 0.988 0.967
Watanabe-Akaike information criterion (WAIC) ...: 36417.72
Effective number of parameters .................: 6.62
Marginal log-Likelihood: -17441.88
is computed
Posterior summaries for the linear predictor and the fitted values are computed
(Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)') mean 0.025quant 0.975quant
lang_fr 1.0857921 1.0591793 1.1120089
lang_it 1.1450520 1.0989781 1.1914942
sexFemale:age_group40-59 0.8912945 0.8495615 0.9350871
sexMale:age_group40-59 1.1642336 1.1221972 1.2078659
sexFemale:age_group60-69 0.9051010 0.8703837 0.9412209
sexMale:age_group60-69 1.2260268 1.1895955 1.2636080
sexFemale:age_group70-79 1.0358802 1.0092955 1.0632110
sexMale:age_group70-79 1.3156365 1.2870635 1.3449175
sexFemale:age_group80+ 1.3131180 1.2948937 1.3318211
sexMale:age_group80+ 0.9852088 0.9699248 1.0008522[1] -17.97158
At first sight, we may thing that there is a large effect of language region on excess mortality, with around 6-11% more deaths than expected in French-speaking municipalities and 10-19% more in Italian-speaking municipalities compared to German. However, as expected this association is likely confounded by the regional variability associated with COVID-19 waves in 2020. Indeed, if we now look at the geographically-structured municipality effect for this model, which can be interpreted as residual effects, we see that the higher excess in South and Southwestern Switzerland is now more evenly distributed (captured by the language effect), while French-speaking regions that were comparatively less impacted during the first wave (such as Neuchâtel and Jura) now have a negative municipality effect to compensate. These nonsensical results highlight the difficulty to estimate the effect of language regions. For this reason, in the following we rely upon observing the residual municipality effects to draw conclusion about the association with language rather than using the language as a fixed effect, as shown in the next map based on model 1.3.
In this last map, we can make two observations. First, French-speaking and Italian-speaking municipalities are not systematically more affected by excess mortality than German-speaking municipalities, with exceptions like the area around Neuchâtel and the Italian-speaking municipalities of Graubunden. Second, there is a clear separation between the French- and German-speaking municipalities, so that differences could rather be attributed to a lower level of connectivity between populations of different language rather than intrinsic differences favoring SARS-CoV-2 transmission or severity.
Referendums on COVID-19 measures
We now focus on results from two referendums about COVID-19 control measures held in June and November 2020. The point here is not to look at causality one way or the other, as we look at overall excess for 2020, and the voting took place at two separated points. A preliminary analysis has reported a negative association between the proportion of “yes” vote at the November referendum at the cantonal level and 7-day incidence on December 7, 2021 (https://smw.ch/index.php/smw/announcement/view/50). We classify municipalities according to the proportion of “yes” vote (expressing support of government-issued measures aimed at controlling COVID-19) at each vote, in quintiles (taking the 5th quintile - with highest support - as a reference).
model1.4e = INLA::inla(munici_observed ~ - 1 + offset(E) +
sex:age_group +
f(id_space, model = "bym2", graph = "data/nb/gg_wm_q.adj", scale.model = TRUE,
hyper = hyper.bym2, constr=TRUE) +
vote_jun_q1 + vote_jun_q2 + vote_jun_q3 + vote_jun_q4,
data = data1,
family = "Poisson",
control.compute = list(config = TRUE, waic = TRUE),
quantiles = c(0.025, 0.5, 0.975),
num.threads = threads,
safe = TRUE)
summary(model1.4e)Time used:
Pre = 21.2, Running = 5.36, Post = 0.48, Total = 27.1
Fixed effects:
mean sd 0.025quant 0.5quant 0.975quant mode kld
vote_jun_q1 0.031 0.015 0.002 0.031 0.061 0.031 0
vote_jun_q2 0.024 0.014 -0.003 0.024 0.051 0.024 0
vote_jun_q3 0.017 0.012 -0.007 0.017 0.041 0.017 0
vote_jun_q4 0.009 0.011 -0.013 0.009 0.030 0.009 0
sexFemale:age_group40-59 -0.095 0.025 -0.144 -0.095 -0.046 -0.095 0
sexMale:age_group40-59 0.172 0.019 0.134 0.172 0.210 0.172 0
sexFemale:age_group60-69 -0.079 0.020 -0.119 -0.079 -0.038 -0.079 0
sexMale:age_group60-69 0.223 0.016 0.191 0.223 0.255 0.223 0
sexFemale:age_group70-79 0.057 0.014 0.029 0.057 0.084 0.057 0
sexMale:age_group70-79 0.295 0.012 0.271 0.295 0.318 0.295 0
sexFemale:age_group80+ 0.294 0.008 0.278 0.294 0.311 0.294 0
sexMale:age_group80+ 0.006 0.009 -0.011 0.006 0.024 0.006 0
Random effects:
Name Model
id_space BYM2 model
Model hyperparameters:
mean sd 0.025quant 0.5quant 0.975quant mode
Precision for id_space 1023.386 243.424 634.025 993.434 1586.815 933.75
Phi for id_space 0.956 0.042 0.842 0.969 0.996 0.99
Watanabe-Akaike information criterion (WAIC) ...: 36434.23
Effective number of parameters .................: 10.42
Marginal log-Likelihood: -17483.54
is computed
Posterior summaries for the linear predictor and the fitted values are computed
(Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)') mean 0.025quant 0.975quant
vote_jun_q1 1.0319485 1.0020940 1.0626900
vote_jun_q2 1.0242462 0.9969791 1.0522754
vote_jun_q3 1.0172793 0.9931611 1.0419920
vote_jun_q4 1.0086063 0.9874339 1.0302451
sexFemale:age_group40-59 0.9092882 0.8659131 0.9548398
sexMale:age_group40-59 1.1874566 1.1432670 1.2333611
sexFemale:age_group60-69 0.9244034 0.8880099 0.9622936
sexMale:age_group60-69 1.2499158 1.2109556 1.2901414
sexFemale:age_group70-79 1.0582487 1.0296109 1.0876954
sexMale:age_group70-79 1.3425417 1.3108470 1.3750249
sexFemale:age_group80+ 1.3423125 1.3205557 1.3644803
sexMale:age_group80+ 1.0064716 0.9886000 1.0246944[1] -1.467352
model1.4f = INLA::inla(munici_observed ~ - 1 + offset(E) +
sex:age_group +
f(id_space, model = "bym2", graph = "data/nb/gg_wm_q.adj", scale.model = TRUE,
hyper = hyper.bym2, constr=TRUE) +
vote_nov_q1 + vote_nov_q2 + vote_nov_q3 + vote_nov_q4,
data = data1,
family = "Poisson",
control.compute = list(config = TRUE, waic = TRUE),
quantiles = c(0.025, 0.5, 0.975),
num.threads = threads,
safe = TRUE)
summary(model1.4f)Time used:
Pre = 22.7, Running = 4.94, Post = 0.489, Total = 28.1
Fixed effects:
mean sd 0.025quant 0.5quant 0.975quant mode kld
vote_nov_q1 0.035 0.015 0.006 0.035 0.064 0.035 0
vote_nov_q2 0.035 0.013 0.009 0.035 0.060 0.035 0
vote_nov_q3 0.025 0.011 0.002 0.025 0.047 0.025 0
vote_nov_q4 0.014 0.011 -0.007 0.014 0.035 0.014 0
sexFemale:age_group40-59 -0.101 0.025 -0.149 -0.101 -0.052 -0.101 0
sexMale:age_group40-59 0.167 0.019 0.129 0.167 0.205 0.167 0
sexFemale:age_group60-69 -0.084 0.020 -0.124 -0.084 -0.044 -0.084 0
sexMale:age_group60-69 0.218 0.016 0.187 0.218 0.250 0.218 0
sexFemale:age_group70-79 0.052 0.014 0.024 0.052 0.079 0.052 0
sexMale:age_group70-79 0.289 0.012 0.266 0.289 0.313 0.289 0
sexFemale:age_group80+ 0.290 0.008 0.273 0.290 0.306 0.290 0
sexMale:age_group80+ 0.002 0.009 -0.016 0.002 0.019 0.002 0
Random effects:
Name Model
id_space BYM2 model
Model hyperparameters:
mean sd 0.025quant 0.5quant 0.975quant mode
Precision for id_space 1115.427 278.225 675.62 1079.557 1763.943 1008.461
Phi for id_space 0.953 0.045 0.83 0.967 0.996 0.989
Watanabe-Akaike information criterion (WAIC) ...: 36432.56
Effective number of parameters .................: 10.09
Marginal log-Likelihood: -17481.17
is computed
Posterior summaries for the linear predictor and the fitted values are computed
(Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)') mean 0.025quant 0.975quant
vote_nov_q1 1.0357006 1.0060449 1.0661876
vote_nov_q2 1.0354068 1.0091542 1.0623005
vote_nov_q3 1.0248533 1.0021952 1.0479715
vote_nov_q4 1.0145273 0.9933728 1.0360981
sexFemale:age_group40-59 0.9043618 0.8612618 0.9496257
sexMale:age_group40-59 1.1813935 1.1375159 1.2269772
sexFemale:age_group60-69 0.9197198 0.8835954 0.9573306
sexMale:age_group60-69 1.2436733 1.2051229 1.2834774
sexFemale:age_group70-79 1.0529810 1.0246059 1.0821643
sexMale:age_group70-79 1.3356940 1.3044409 1.3677354
sexFemale:age_group80+ 1.3358507 1.3145561 1.3575778
sexMale:age_group80+ 1.0015842 0.9840605 1.0194700[1] -3.13754
In both cases it appears that excess mortality is consistently about 1-6% higher in municipalities expressing lowest support to control measures (first quantile) in both referendums, with a clear gradient.
Multivariable model
We now jointly estimate the effects of interest identified in the univariable analysis: rural or urban status, border, median SEP quintile and results from the COVID-19 referendums (using only the June referendum to limit complexity). We leave out language regions for the reasons explained above.
model1.5 = INLA::inla(munici_observed ~ - 1 + offset(E) +
sex:age_group +
f(id_space, model = "bym2", graph = "data/nb/gg_wm_q.adj", scale.model = TRUE,
hyper = hyper.bym2, constr=TRUE) +
border +
type_periurban + type_urban +
sep1 + sep2 + sep3 + sep4 +
vote_jun_q1 + vote_jun_q2 + vote_jun_q3 + vote_jun_q4,
data = data1,
family = "Poisson",
control.compute = list(config = TRUE, waic = TRUE),
quantiles = c(0.025, 0.5, 0.975),
num.threads = threads,
safe = TRUE)
summary(model1.5)Time used:
Pre = 20.3, Running = 5.72, Post = 0.659, Total = 26.6
Fixed effects:
mean sd 0.025quant 0.5quant 0.975quant mode kld
border 0.028 0.013 0.002 0.028 0.054 0.028 0
type_periurban -0.005 0.013 -0.031 -0.005 0.020 -0.005 0
type_urban -0.019 0.013 -0.044 -0.019 0.006 -0.019 0
sep1 0.061 0.017 0.027 0.061 0.095 0.061 0
sep2 0.036 0.014 0.008 0.036 0.064 0.036 0
sep3 0.031 0.013 0.005 0.031 0.056 0.031 0
sep4 0.003 0.011 -0.018 0.003 0.025 0.003 0
vote_jun_q1 -0.008 0.018 -0.043 -0.008 0.027 -0.008 0
vote_jun_q2 -0.006 0.016 -0.037 -0.006 0.025 -0.006 0
vote_jun_q3 -0.005 0.013 -0.031 -0.005 0.021 -0.005 0
vote_jun_q4 -0.004 0.011 -0.026 -0.004 0.017 -0.004 0
sexFemale:age_group40-59 -0.095 0.029 -0.151 -0.095 -0.039 -0.095 0
sexMale:age_group40-59 0.172 0.024 0.125 0.172 0.219 0.172 0
sexFemale:age_group60-69 -0.078 0.025 -0.126 -0.078 -0.029 -0.078 0
sexMale:age_group60-69 0.223 0.021 0.181 0.223 0.265 0.223 0
sexFemale:age_group70-79 0.057 0.020 0.018 0.057 0.095 0.057 0
sexMale:age_group70-79 0.295 0.018 0.259 0.295 0.331 0.295 0
sexFemale:age_group80+ 0.295 0.016 0.263 0.295 0.327 0.295 0
sexMale:age_group80+ 0.007 0.017 -0.026 0.007 0.039 0.007 0
Random effects:
Name Model
id_space BYM2 model
Model hyperparameters:
mean sd 0.025quant 0.5quant 0.975quant mode
Precision for id_space 1423.440 408.983 801.057 1363.295 2397.106 1246.678
Phi for id_space 0.939 0.065 0.758 0.961 0.997 0.992
Watanabe-Akaike information criterion (WAIC) ...: 36429.45
Effective number of parameters .................: 10.52
Marginal log-Likelihood: -17526.34
is computed
Posterior summaries for the linear predictor and the fitted values are computed
(Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)') mean 0.025quant 0.975quant
border 1.0283387 1.0016133 1.0558011
type_periurban 0.9945418 0.9690752 1.0206705
type_urban 0.9812836 0.9568599 1.0063273
sep1 1.0627038 1.0270331 1.0993785
sep2 1.0367913 1.0084486 1.0657314
sep3 1.0310638 1.0051234 1.0575161
sep4 1.0033027 0.9817893 1.0252365
vote_jun_q1 0.9918056 0.9576895 1.0272466
vote_jun_q2 0.9938405 0.9635376 1.0251943
vote_jun_q3 0.9951424 0.9694527 1.0215967
vote_jun_q4 0.9955408 0.9741797 1.0174281
sexFemale:age_group40-59 0.9095203 0.8599219 0.9619987
sexMale:age_group40-59 1.1876306 1.1334772 1.2444052
sexFemale:age_group60-69 0.9253048 0.8813274 0.9715011
sexMale:age_group60-69 1.2498827 1.1988881 1.3030943
sexFemale:age_group70-79 1.0582221 1.0180487 1.1000260
sexMale:age_group70-79 1.3432789 1.2954773 1.3929092
sexFemale:age_group80+ 1.3425998 1.3004781 1.3861746
sexMale:age_group80+ 1.0065400 0.9743249 1.0398804[1] -6.244418
This multivariate analysis confirms the association between high relative excess mortality and low median SEP, with a clear gradient. We also see that the association between high excess mortality and border municipalities remains, and than urban areas are still associated with comparatively higher excess mortality in 2020, although the data is now compatible with no effect. Estimates of an association with voting results are faint after adjusting for the other covariates. This decrease in the association can be attributed to the collinearity with other covariates, in particular with SEP and urbanicity (see correlogram in section 1). Therefore, we cannot disentangle between the associations with voting results, SEP and urbanicity.
There are also some interesting patterns in the residual effects at the level of the municipality (adjusting for all aforementioned covariates), with in particular, expected higher excesses in Ticino and Southwestern Switzerland, a visible barrier between French-speaking and German-speaking regions, lower excess in the large cities of the German-speaking part (Zurich, Basel, Bern) and in relatively isolated valleys of Graubunden.
Step 2: multivariate model with full uncertainty propagation
In the previous section we modeled the variation of the median excess mortality over 2020 by municipality. This approach underestimates the uncertainty from two sources, first from the prediction error in the expected mortality at the cantonal level, second from the downscaling to the municipal level. At this stage, we bring back these two sources of uncertainty in the final estimates by repeatedly fitting model 1.5 to 50 different sets of posterior draws of excess mortality by municipality, then combining with equal weights.
Time used:
Pre = 953, Running = 210, Post = 22.3, Total = 1185
Fixed effects:
mean sd
border 0.031 0.023
type_periurban 0.004 0.021
type_urban 0.008 0.020
sep1 0.043 0.027
sep2 0.030 0.023
sep3 0.023 0.022
sep4 0.005 0.018
vote_jun_q1 -0.015 0.031
vote_jun_q2 -0.014 0.029
vote_jun_q3 -0.008 0.023
vote_jun_q4 -0.010 0.019
sexFemale:age_group40-59 -0.186 0.052
sexMale:age_group40-59 0.137 0.043
sexFemale:age_group60-69 -0.149 0.051
sexMale:age_group60-69 0.179 0.043
sexFemale:age_group70-79 0.029 0.042
sexMale:age_group70-79 0.248 0.038
sexFemale:age_group80+ 0.277 0.036
sexMale:age_group80+ -0.011 0.039
Random effects:
Name Model
id_space BYM2 model
Model hyperparameters:
mean sd
Precision for id_space 170.283 56.882
Phi for id_space 0.616 0.201
Marginal log-Likelihood: -41448.59
is computed
Posterior summaries for the linear predictor and the fitted values are computed
(Posterior marginals needs also 'control.compute=list(return.marginals.predictor=TRUE)')
As expected, this approach leads to a dilution of the observed associations between relative excess mortality and local covariates. We still observe a linear gradient in the association between excess mortality and median SEP at the municipal level, and a likely association between excess mortality and border municipalities, although in both cases with higher uncertainty. We also observe similar patterns in the residual effects at the level of the municipality, again with higher uncertainty.